﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;

namespace Nature.Toolkit.GKin
{
    public abstract class SoneStatePostUpdateChemReactorController :
        StatePostUpdateChemReactorController, 
        Nature.Numerics.ISoneFunction,
        Nature.Numerics.ISoneVariableCorrection
    {
        private Nature.Numerics.Basic.SoneNewton<SoneStatePostUpdateChemReactorController> m_newton;
        private double[] m_yk;
        private object m_anchor;
        private const double Ymin = 1.0e-30;

        protected override void Reset()
        {
            if (m_newton == null || m_newton.Nva != this.Nsp)
            {
                m_newton = new Numerics.Basic.SoneNewton<SoneStatePostUpdateChemReactorController>(this);
                m_yk = Realloc(m_yk, Nsp);
                m_newton.SetAbsoluteConvergenceRadius(1.0e-4);
                m_newton.SetRelativeConvergenceRadius(1.0e-1);
                m_newton.SetAbsoluteTolerance(1.0e-8);
                m_newton.SetRelativeTolerance(1.0e-5);
            }
        }

        protected override double AdvanceChemicalComposition()
        {
            var yk = MixtureSavePoint.SpeciesMassFractions;
            for (int isp = 0; isp < Nsp; ++isp)
            {
                m_yk[isp] = yk[isp];
            }
            m_newton.Solve(m_yk, ref m_anchor);

            return base.NextTimeIncrement;
        }

        protected override void Rollback()
        {
            base.Rollback();
            ProposeNextTimeIncrement(0.5d * base.NextTimeIncrement);
        }

        protected override void Commit()
        {
            base.Commit();
            ProposeNextTimeIncrement(1.1d * base.NextTimeIncrement);
        }

        #region ISoneFunction Members

        int Numerics.ISoneFunction.Nva
        {
            get { return Nsp; }
        }

        void Numerics.ISoneFunction.Compute(double[] speciesMassFractions, double[] equation)
        {
            var m0 = base.MixtureSavePoint;
            var m1 = base.Mixture;
            m1.UpdateTPY(m0.Temperature, m0.Pressure, speciesMassFractions);
            var yk0 = m0.SpeciesMassFractions;
            var yk = m1.SpeciesMassFractions;
            var wk = m1.SpeciesMassFractionChemicalNetProductionRates;
            var dt = this.NextTimeIncrement;
            for (int isp = 0; isp < Nsp; ++isp)
            {
                equation[isp] = (yk[isp] - yk0[isp]) - dt * wk[isp];
            }
            
        }

        #endregion

        #region ISoneVectorCorrection Members

        void Numerics.ISoneVariableCorrection.ApplyCorrection(double[] x)
        {
            double tmp;
            double sum = 0d;
            for (int isp = 0; isp < Nsp; ++isp)
            {
                tmp = x[isp];
                if (tmp < 0d) { tmp = Ymin; x[isp] = tmp; }
                sum += tmp;
            }
            for (int isp = 0; isp < Nsp; ++isp)
            {
                x[isp] /= sum;
            }
        }

        #endregion
    }
}
